Matlab
code 3.1: Matlab file “Figure 3-6.m”
%--------------------------------------------------------------------
% This
code can be used to generate Figure 3.6
%--------------------------------------------------------------------
clear
all;
close
all;
%% the
theoretical curve of the rotating target micro-Doppler characteristic(LFM
radar)
c = 3e8;
j =
sqrt(-1);
fc =
10e9; % carrier frequency of transmitted signal
B =
200e6; % bandwidth
cord =
1000*[3 4 5]; % coordinates of local coordinate system's origin in the radar
coordinate system
colo =
[0 0 0]; % coordinates of radar in the radar coordinate system
R0 =
cord-colo;
n =
R0/sqrt(sum(R0.^2)); % unit vector of Line-of-Sight(LOS)
ae = [0
pi/4 pi/5]; % Euler angle
tgt = [0
0 0;3 1.5 1.5;-3 -1.5 -1.5]; % three point scatterers
w = [pi
2*pi pi]; % angular velocity of rotation
T =
0.8165; % rotating period
ri1 =
[cos(ae(3)) sin(ae(3)) 0;-sin(ae(3)) cos(ae(3)) 0;0 0 1];
ri2 = [1
0 0;0 cos(ae(2)) sin(ae(2));0 -sin(ae(2)) cos(ae(2))];
ri3 =
[cos(ae(1)) sin(ae(1)) 0;-sin(ae(1)) cos(ae(1)) 0;0 0 1];
ri =
ri1*ri2*ri3; % initial rotation matrix
tn =
length(w); % number of scatterers
r0 =
zeros(3);
for i =
1:tn
r0(i,:) = tgt(i,:)-colo; % scatterers in
the local coordinate system
end
r0r =
ri*r0'; % scatterers in the reference coordinate system
w =
ri*w'; % angular velocity of rotation in the reference coordinate system
omega =
sqrt(sum(w.^2));
we =
w/omega; % unit vector of rotation
wr = [0
-we(3) we(2);we(3) 0 -we(1);-we(2) we(1) 0]; % skew symmetric matrix
v = 0; %
translational velocity of target
t = 3; %
radar illumimated time
tp =
5e-6; % pulse duration
prf =
1000; % pulse repetition frequency
pri =
1/prf; % pulse repetition interval
dt =
0:pri:t-pri; % slow time sampling interval
m =
length(dt);
rdt =
zeros(tn,m);
for i =
1:m
rdt(:,i) =
((eye(3)+wr*sin(omega*dt(i))+wr^2*(1-cos(omega*dt(i))))*r0r)'*n';
end
plot(dt,rdt)
xlabel('Slow
time (s)')
ylabel('Range
(m)')
axis([0,3,-6,6])
%% the
simulation result of the rotating target micro-Doppler characteristic(LFM
radar)
xmin =
sqrt(sum((cord-colo).^2)); % the minmum detectable distance
xmax =
sqrt(sum((cord-colo).^2))+sqrt(sum((tgt(2,:)-colo).^2)); % the maximum detectable
distance
ts =
xmin/c; % the starting point of the time sampling
te =
xmax/c+tp; % the ending point of the time sampling
mu =
B/tp; % frequency modulation rate
fs = B;
% fast time sampling frequency
kdt =
1/fs; % fast time sampling interval
nt =
2*ceil((te-ts)/(2*kdt)); % number of fast time sampling
df =
fs/nt; % number of frequency sampling
Rr =
zeros(length(tgt(1,:)),length(dt)); % distance between the scatterers and radar
s =
zeros(nt,length(dt)); % echo signal matrix
for i =
1:m
for n = 1:length(tgt(:,1))
Rr(n,i) =
sqrt((R0'+v*dt(i)+((eye(3)+wr*sin(omega*dt(i))+wr^2*(1-cos(omega*dt(i))))*r0r(:,n)))'...
*(R0'+v*dt(i)+((eye(3)+wr*sin(omega*dt(i))+wr^2*(1-cos(omega*dt(i))))*r0r(:,n))));
end
end
tk = ts+kdt*(0:nt-1);
% fast time
for n =
1:length(tgt(:,1))
td =
tk(:)*ones(1,m)-2*ones(nt,1)*Rr(n,:)/c; % echo delay
s = s+exp(j*2*pi*fc*td+j*pi*mu*td.^2);
end
Rr0 =
sqrt(sum(R0.^2));
td0 =
tk(:)*ones(1,m)-2*Rr0/c;
s0 =
exp(j*2*pi*fc*td0+j*pi*mu*td0.^2); % reference signal
sdt =
s.*conj(s0); % dechirp
sdf =
fftshift(fft(sdt),1);
x = dt;
y =
-(df*(-nt/2:nt/2-1))*c/(2*mu); % turn the fast time frequency-slow time plane
to the range-slow time plane
figure
contour(x,y,abs(sdf),10)
xlabel('Slow
time (s)')
ylabel('Range
(m)')
axis([0,3,-6,6])